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. stochastic Langevin type heat bath at each site. The temperatures 

I of the "exterior" left and right heat baths are at specified values 

■ and Tr, respectively, while the temperatures of the "interior" baths 

are chosen self-consistently so that there is no average flux of energy 
between them and the system in the steady state. We prove that 
this requirement uniquely fixes the temperatures and the self consis- 
tent system has a unique steady state. For the infinite system this 
state is one of local thermal equilibrium. The corresponding heat cur- 
rent satisfies Fourier's law with a finite positive thermal conductivity 
which can also be computed using the Green-Kubo formula. For the 
^ , harmonic chain (d = 1) the conductivity agrees with the expression 

^ [ obtained by Bolsterli, Rich and Visscher in 1970 who first studied 

this model. In the other limit, d ^ 1, the stationary infinite vol- 
ume heat conductivity behaves as {£dd)~^ where id is the coupling 
to the intermediate reservoirs. We also analyze the effect of having 
a non- uniform distribution of the heat bath couplings. These results 
are proven rigorously by controlling the behavior of the correlations 
in the thermodynamic limit. 
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1 Introduction 

Our understanding of non-equilibrium systems is at the present time very 
incomplete. In particular, we still have no model Hamiltonian system for 
which Fourier's law has been proven rigorously. A review of the problems 
and of the few known exact results related to Fourier's law is given in and 
in [2] which also contains a survey of recent numerical results. 

Here we study, in a mathematically rigorous manner, the microscopic 
structure of the stationary non-equilibrium state of a "self-consistent har- 
monic crystal", a model introduced by Bolsterli, Rich and Visscher (BRV) in 
[3 0] . This is a d-dimensional system of iVi x ■ • • x A^^ oscillators whose time- 
evolution is given by a combination of Hamiltonian and stochastic dynamics. 
The Hamiltonian is composed of harmonic nearest neighbor and "on-site" 
potentials and the stochastic part comes from coupling each particle in the 
chain to its own heat bath. 

We want to describe a situation where we have a temperature gradient 
in one direction (the "first" with A^i oscillators) while the temperature is 
uniform in the remaining d — 1 directions, on which we impose periodic 
boundary conditions. The temperatures of the heat baths of the end-point 
particles in the first direction are fixed to given values Tl and Tr, while 
the temperatures of the interior heat baths are chosen self-consistently by 
the requirement that there is no energy flux, on average, between any such 
reservoir and the system in the steady state. From a physical point of view, 
we may think of the interior heat reservoirs as representing schematically the 
effect of degrees of freedom not included in the Hamiltonian. 

Using numerical studies and non-rigorous arguments, BRV found that in 
the case d = 1 [chain) the (kinetic) temperature profile of the system in 
its steady state is linear, with a heat flux proportional to for large N. 
This corresponds to the self-consistent system having a finite, temperature 
independent, thermal conductivity P3I2]- These results are in sharp contrast 
to those found earlier by Rieder, Lebowitz and Lieb 0, who studied a system 
with the same Hamiltonian dynamics, but with heat baths acting only at the 
boundaries. They found that the system had an infinite conductivity and a 
constant temperature profile away from the ends, results later generalized to 
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the higher dimensional case by Nakazawa [Hj- 

In this paper, we provide a rigorous proof that the steady state of the 
self-consistent system has indeed the properties found by BRV for the d = 1 
case in 13 H], and we extend the results to cover all d > 1. More precisely, 
we show that, in the limit where all Ni —>■ oo, the steady state is a local 
equilibrium state j7] with a temperature profile satisfying Fourier's law with 
a finite, temperature independent, thermal conductivity. 

We deal first with the d = 1 case, and consider the higher dimensional 
case only in section [7| We define the model and solve its dynamics with a 
given temperature profile in section |21 We then turn to the self-consistency 
condition in section El proving in particular, that the self-consistent profile 
is always uniquely determined by the boundary temperatures. Section HI 
contains our main results: we prove there the local equilibrium property and 
Fourier's law. In section we use the explicit solution to show that the 
Green-Kubo formula holds for this system. In section IHl we briefly analyze 
the case where the couplings to the heat baths are non-uniform, and we 
conclude that the local macroscopic heat conductivity is proportional to the 
inverse of the local average of the couplings. 

Finally in section[7J we first show how one can map the higher dimensional 
self-consistent system with periodic boundary condition on all directions but 
the first into a set of one- dimensional chains, and then apply the earlier results 
to derive a generalization to the higher dimensional case. We show there that 
for large d the conductivity behaves as {idd)~^ where £d is the coupling to 
the intermediate reservoirs. Some technical details of the calculations are 
collected in the appendices. 

2 Dynamics and the stationary state 

To start with, we consider a chain of N oscillators with a Hamiltonian 



where q and p are vectors in R , we set go = = Qn+i, and we have 



N , N+1 



(2.1) 



i=l i=l 



u{q) = -7^g^ and v{q) = 




(2.2) 
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with 7, > 0. In addition, the oscillator at each site i is coupled to a 
Langevin heat bath at temperature Tj > 0, with a coupling strength A > 0. 
As in the time-evolution of the system is then given by the stochastic 
differential equations, 

X = -AX + J:W (2.3) 

where X = {q,p) is the phase-space vector, the Wi are independent white 
noises and A and E are 2N by 2N matrices given by 

aJI -n a„d E^(° ° V (2.4) 



$ A y V v/2Ar 

Here / is the unit by matrix, Aij = SijX, %j = SijTi, and 

^ = uj\-A + u^I), (2.5) 

where z/^ = and A denotes the discrete Laplacian with Dirichlet bound- 

ary conditions: 

Ajj = —26ij + + 

Equations ()2.3|) define an Ornstein-Uhlenbeck process whose solution with 
initial data X(0) is given by the stochastic integral (for details, see e.g. 
chapter 5 in [8j) 

X{t) = e''^X{0) + rdse-(*-^)^EVr(s). (2.6) 
Jo 

This is a Gaussian process, determined uniquely by its mean and covari- 
ance which can be computed directly from (j2.6p . However, for latter use we 
assume now that X(0) is distributed according to a Gaussian measure with 
a mean Xq and a covariance Co — the deterministic case is then obtained by 
setting Co = 0. Then the mean evolves by 

(X(t)) = e-*^Xo, (2.7) 

and for the covariance C{t',t) = {[X{t')-{X{t'))] ® [X{t)-{X{t))]) we get 

^^^^^>-\ C{t',t')e-('~'')^^ iit'<t ' ^^-^^ 
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with 

C(t,t) =e-*^Coe-*^'^+ Tdse-^^EV^^^. (2.9) 

Jo 

We show in Appendix El that for any a satisfying 



A 7^ 

< a < min<i — , — 



2' A 

we can find a constant c < oo such that for alH > and for all N 

||e~*^|| < ce"*". (2.10) 

The uniform exponential decay of e"*"^ implies that there is an exponentially 
fast convergence in the microscopic scale to a unique stationary state, which 
is Gaussian with mean and covariance S, 

POO 

S= / dse-^^SV^^''. (2.11) 
Jo 

This S is the unique solution of the equation 

AS + 3X^ = 1:^, (2.12) 

which we solve following ref. 0. We divide S into by components, 



S 



u z 

V 



and get the following four equations equivalent to (j2.12j) : 

Z = -Z^, r = ^($f/ + [/$) + ^(ZA-AZ), 

ZK^ KZ = ^u -u^, ^^■^'^> 
A(r - K) + (T - V)k = ^Z - z^. 

A diagonalization of the discrete Laplacian yields 



5 



where Y^i = jj-k^kh ^-k are the eigenvalues of 



^ = z/2 + 4sin2f— (2.14) 



and F is the orthonormal matrix 



, 2 / Txkl , 

^'i' = \/l^^^sin — — . (2.15) 



A^+1 VA^ + 1 

As shown in Appendix ^ any block B of the covariance matrix (i.e. B = U, 
V or Z) can then be obtained by a linear transformation of the form 

N 

B^J=Y.4^^r (2.16) 



r=l 



where 



AT 



^ij^ = Yl F^kF,lf^''\ck, Cl)FrkFrl (2.17) 



k,l=l 



Cfc = COS ( — — ) . (2.18) 



and 

r.i^ = nos I 

^A^+ 1 

The functions f^^^ for the different choices of B are given by 

A2 1 



uj^Gix,y) 



where 



G{x, y) = (x-yf + ^iu^ + 2-x-y). (2.20) 

When Ti = T for all i, only the values with k = I contribute in equation 
(|2.17p . The stationary covariance is then given by 

= T ' (2.21) 

and the stationary measure is the Gibbs measure at temperature T. 
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2.1 Energy current 

We define the local energy of particle i by 

Hi{q,p) = ]^p1 + u{qi) + ^ [v{qi - qi_i) + v{qi+i - g^)] 

for z = 2, . . . , — 1. The boundary terms {i = 1, A^) are defined similarly, 
but with double the interaction energy contribution from the connections to 
go = and to ^at+i = 0. With these definitions, 

N 

H{q,p) = Y,H,{q,p) 

i=l 

and 

^^{HM = - mt)) - {J^-lm + (2.22) 



where 



Ji{QiP) = — ^{Qi+i - Qt) {Pi + Pi+i) , i = l,...,iV- 1, (2.23) 
Jj = for i = 0, A^, and 

R,{q,p) = X{T,-p'^). 

The Ji correspond to the energy currents inside the system while Ri gives 
the energy flux from the i-th reservoir to the i-th oscillator. 

The corresponding expectation values in the stationary state are 

^•^'^^ = Y - + + = ^'^..+1, (2-24) 

where we have used the antisymmetry of Z, and 

{Ri)s = \{Ti- Vii) . (2.25) 



3 Self-consistency condition 

As described in the introduction, we let the end-point temperatures Ti = Tl 
and T/v = Tr be given independently of A^. Then we want to choose Ti for 
i = 2, . . . , A^ — 1 in such a way that {Ri)s = or, by ()2.25p . so that 

T, = V,i. (3.1) 
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By (|2.1ip . the kinetic temperature vector (Vij)^]^ depends linearly on 
the imposed temperature vector T. Therefore, all solutions to the self- 
consistency condition can be obtained using the equation 

t, = Tr + {Tl - rR)rf (3.2) 

where T^^'°^ is a solution of the problem with Tl = 1 and Tr = 0. The set of 
all such T^^'^^ form a convex set, i.e. any non-uniqueness in the solution of the 
self-consistency condition would imply the existence of a whole continuum 
of solutions. We shall now prove that for our model the solution to the 
self-consistency problem is unique. 

Let us denote the above linear mapping from imposed to kinetic tempera- 
tures by M, i.e. Vu = Ylj ^ij^j- It follows straightforwardly from ()2.11|) that 
Mij > for all pairs i,j, and by the explicit solution given in ()2.16|) - ()2.20p 
we have 

N 

^ij = '^iP = FikFiifkiFjkFji, (3.3) 

k,l=l 

where 

/.z = /(^)(c.,Q) = l-^?^^>0. 

G{ck,ci) 

By ()3.3|1 . M is then symmetric and satisfies for all i 

N N 

5^M,, = 5^M,, = 1, (3.4) 

which imply that M is, in fact, a doubly stochastic matrix. 

Theorem 3.1 For any given end-point temperatures Tl, > 0, there is a 
unique, positive temperature profile which satisfies the self-consistency con- 
dition. In addition, all temperatures in the profile lie between the end-point 
temperatures. 

Proof: Let T be a self-consistent profile and define 

T,-Tn, forz = 2,...,Ar-l, 
0, if z = 1 or i = 

Clearly a belongs to the subspace 

W = IxeR^ xi = = XN 
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Let us denote the orthogonal projection to the subspace W by Pw, and define 
Q = PwMPw. Using (j3.4|) in (j3.H) we get the equation 

a = Qa + b (3.5) 

where the vector b is defined by 



b. 



MaiTL-TR), for^ = 2,...,iV-l, 
0, if i = 1 or i = 



We shall later prove in Corollarv 13 . 31 that \\Q\\ < 1. Then for any b G M^, 
equation ()3.5|) has a unique solution given by 

oo 

a = Y, Q" b, 

n=0 

and the self-consistent profile must thus satisfy 

oo 

Ti=TR + Y,{Q''b)u for ^ = 2, . . . , iV - 1. (3.6) 

n=0 

As the profile defined by (j3.6|) also satisfies (|3.1|) . this proves the existence 
and uniqueness of the self-consistent profile. Repeating the above computa- 
tion for a' = PwiT — T^l) instead of a, and then using the positivity of Mjj 
in the equations corresponding to ()3.6p . proves that Tr < Ti < Tl for all i 
when Tr <Tl, and that <Ti < Tr when Tr>Tl. □ 

To conclude the proof, we still need to prove that \\Q\\ < 1. This will 
follow from the following lemma: 

Lemma 3.2 Let c he a constant which satisfies 

c 



for all X, ?/ G [— 1, 1]. Then, for any vector x G 



pN 



x^(I - M)x>-\\Dxf, (3.7) 

c 



where D is the 'finite difference operator" defined by 
{Dx)i = 



Xi — Xj+i, for 1 < i < N — 1 
0, fori = N 
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Proof: Let x G M^, and define Xki = J2i=i FikFuXi for all fc, /. Then by 
(13. 3|) and the assumption G{ck,ci) < c/2, 

2 ^ 

x^{I-M)x>-J2S:Uck-ciy 
c ^ — ^ 

k,l=l 
2 JV TV 

= - ^ XiXj ^ FikFiiFjkFji{{l - Ck) - (1 - q))^. 

i,j=l fc,/=l 

But by an explicit computation, we have for all i, j, 

N 

2 FikFuF,kFji{{l - Ck) - (1 - cz))2 = (A^),,.^,,. - {Ai,f = {D^D)i„ 

k,l=l 

and the inequality ()3.7|) has been derived. □ 

For the applications, we note that a constant c satisfying the requirement of 
the lemma can always be found. 

Corollary 3.3 ||Q|| < 1. 

Proof: Since W is finite-dimensional and Q = PwQPw is symmetric, it 
will be enough to show that for ail x eW with ||a;|| < 1 

< x^Qx < 1. (3.8) 

But for such x, x^Qx = x'^Mx, and the lemma together with fki > yields 

< x'^Mx < 1 - x^{I - M)x < 1 - ill/^icll^ (3.9) 

If Dx = 0, then = for all i, and we must now have either > 0, 

when fj3.9|) implies ()3.8|) . or a:; = 0, when ()3.8|) is trivially true. □ 

4 Fourier's law 

In this section we first derive a number of technical estimates which will be 
necessary to control the behavior of the system in the thermodynamic limit. 
We then show that the self-consistent steady state is microscopically a local 
equilibrium state with a heat flux satisfying Fourier's law. 
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4.1 Exponential decay of correlations 

Let / denote any one of the three functions f^^^ defined by equation ()2.19p . 
Clearly, / is a rational function, analytic everywhere but at the zeroes of G 
in flT^ . On the other hand, 

G{x,y)>^ = X'^^ (4.1) 

for all X, y E [—1, 1], and since we have assumed that 7 > 0, there are no 
zeroes of G inside [—1,1]^. 

This implies that the function /(cos ■, cos ■), which enters in ()2.17j) . is 
analytic in some neighborhood region of [— 7r,7r]^ in C^. In particular, its 
Fourier series converges pointwise, and we have for all x, ?/ G M, 

00 

/(cosx, cosy) = cos(mx) cos(ny)/(m, n). (4.2) 

m,n=~oo 

Applying this with / = f^^^ in ()2.17p we get after some straightforward 
algebra, 

B^j^ = fN{i-r,j-r) + fN{i + r,j + r)~fN{i-r,j + r)-fN{i + r,j~r) (4.3) 

where is defined by 

Mm, n)=J2 /("^ + 2(A^+1)A;, n + 2{N+l)l). (4.4) 
k,lez 

By the above mentioned analyticity, the Fourier coefficients / decay ex- 
ponentially. From this the following behavior for is easily derived: 

Lemma 4.1 For any block B, define J'n by ^.4\ ) using f = f^^^. Then 
there are strictly positive, N -independent constants a and a such that, for 
any m,n & Z, 

|/7v(m,n)| <ae-'^(l'"'l+l"'l) (4.5) 
where m' = (m mod 2(A^+ 1)) G {—N, . . . , 1} and n defines n' similarly. 



11 



4.2 Local equilibrium 

Consider some temperature profile T = (Tj)^]^, and let its maximum nearest 
neighbor variation be denoted by en, i-e. with D defined as in Lemma l3.2[ 
let 

En = max \ {DT)i\ = max |Tj — Tj+i|. (4.6) 

i i<N 

We shall prove in this section that, if <^ 1 and maxj Tj is bounded, the 
local microscopic properties of the stationary measure near a site i can be 
well approximated by using the equilibrium measure with the temperature 
Tj. This will, in particular, justify our identification of the parameter Tj as 
a local temperature. 

The main ingredient of the proof is the following corollary to Lemma [4. II 



Corollary 4.2 For any block B, there is an N -independent constant a' , such 
that for all i,j G {1, . . . , A^}, 



N 



^|z-r||i?J)|<a'. (4.7) 



r=l 



Proof: Apply equation (j4.3|) . the triangle inequality, and Lemma f4. II to re- 



place \bIj^\ in (j4.7|) by four exponential terms. The estimate follows straight- 
forwardly. □ 

Let B be any block of the stationary covariance matrix, and let i, j & 
{1, . . . , A^}. As it takes |r — 2| "steps" to get from a site r to the site i, we 
have the obvious bound 

\Tr-Ti\ <\r -i\EN. (4.8) 
Then Corollarv 14.21 immediatelv yields the estimate 

N 

B,,-T,J2Mj^\ <o^'^N. (4.9) 

r=l 

The value of Xl^Li -^i]'' gives the component of the stationary covari- 
ance matrix when all the temperatures are set equal to Tj, i.e. the equilibrium 
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covariance at temperature T^, see ^TT^ . As a' in (j4.9p can be chosen 

independently of A^, we have now shown that all pair-correlations satisfy 

B,,=BlP'^^^ + 0{eN). (4.10) 

As both the equilibrium measure and the above stationary measure are 
Gaussian, we can also conclude that all finite correlations can be approxi- 
mated by the local equilibrium values, with an error which vanishes when 
En 0. 



4.3 Self-consistent current and Fourier's law 

Let T now denote the self-consistent profile which by Theorem 13.11 is unique 
and is bounded by Tl and Tr. Then by fl2.13p and ()2.24|) we get for all 
i = 2,...,N-l, 

= A(r - Vh = \{^Z - Z^\, = a;2(-Z,„i,, - Z,+i,,) = {J,)s - {J^-l)s. 
Therefore, the steady state current is constant throughout the chain: 

(J,)s = /^\ for^ = l,...,Ar-l. (4.11) 

The same equations also imply the following relation between the local 
steady state current and the q-q correlations: for alH = 1, . . . , A^ — 1, 

{J^)s = uj^Z,,+r = ^($ _ U<^\,+i 

(4.12) 



where we define Uij = when i or j ^ {1, . . . , A^}. Summing ()4.12|) over all 
the indices i, we get by using equation (j4.11|) 

(Ar-l)jW = |^(f/n-^^;v). (4.13) 

Let then D and c be given as in Lemma 13. 2^ and as before let = 
maxj \ {DT)i\. Then the magnitude of the current and En are related by the 
formula 



eN<\ljiTL~Tn)J(^l (4.14) 
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To see this, first let x = {I — M)T where the matrix M was defined in 
sectional Since T is self-consistent, Xj = except possibly at the end-points. 
But then by ()3.4|) . xi+xn = Yli = O5 thus also xn = —xi. Therefore, 
fj4.14p follows from Lemma f3. 21 since then 



xi{Tl - Tr) = T^{I - M)T > -\\DTf > -e 



2 



and, by x^=T^- V^n = U^Z - Z<l>)n = (Ji)/A. 

Using ()4.14j) we can estimate the error made when the terms on the right 
side of fj4.13|) are replaced by their equilibrium values. Equation ()4.9|) and 



the explicit form of the equilibrium covariance given in ()2.2H1 yield 



,2 



<a'yJ^(Tz.-r^)jW (4.15) 



where a' is chosen as in Corollarv 14.21 for the block B = U . 

It follows from symmetry that = '^'n^, which has the limit 

T 2 sin^fTra;) 2 1 , , 

lim $7/ = — / dx— ^ ' , = — , 4.16 

Jo z/2 + 4sin2(f) uj^2 + v^ + ^v\A + v^) ^ ' 

It is then a consequence of ()4.15j) that 

lim (iV - 1) J(^) = k{Tl - Tr) (4.17) 

where 

K = ^ \ (4.18) 

A 2 + z/2 + v/z/2(4 + z/2) 

Since 

JW = O(A^-i), we get from (imil that 

ej^ = 0{N-'2). (4.19) 

By our discussion in section 14.21 this implies that in the limit N ^ 00, 
all the correlation functions involving finitely many terms will converge to 
the corresponding local equilibrium values if we identify Tj with the local 
temperature of the system at a site i. 

Summing ()4.12|1 from 1 to j — 1 and combining it with ()4.13|1 yields 

^^(t/11 - f/iv^) = Uu - U,, + f/,_i,,+i. (4.20) 

14 



Then the local equilibrium approximation, equation (j4.1Up . shows that 

^^<l>r/(T^ - Tn) = n^u' - Tj^J- + T,-i^J-w + 0(6,,). (4.21) 
For all j and k, Fj^i^kFj+i^k = F'^k ~ -^u^ ^^'^ t\ms satisfies the identity 

= ^]}- ^-.l (4.22) 
Since < 1/7^ and > uniformly in A^, we then get the result 

where the correction term vanishes uniformly in j when N —>■ oo. 

Setting X = j/N the system therefore approaches, in the limit N oo, 
a local equilibrium state with a temperature profile 

T{x)=TL + xiTR-TL), xe[0,l]. (4.24) 

Thus Fourier's law holds for the steady state of the system, and the thermal 
conductivity is given by the temperature independent constant k in (j4.18p . 
Note that k, remains finite when z/ — > 0, which points towards a finite con- 
ductivity even for the system without the on-site binding potential. 

Let us finally remark that we do not think the above bound for e^v is op- 
timal. Preliminary numerical simulations suggest that = 0{N^^) rather 
than 0{N-^). 



5 The Green-Kubo formula 

The Green-Kubo formula expresses the local equilibrium conductivity at a 
position X with temperature T{x) as an integral over the current-current 
correlations in a (closed) equilibrium system at uniform temperature T = 
T{x). This corresponds, for the type of stationary state we consider, to a 
formula for k, when and Tr T. It is not immediately apparent how the 
presently available derivations of such a formula (for recent results, see e.g. 
inmni) couM be appHed to a stochastic system like the one considered here. 
In particular, it is not clear which current we should use in the formula: i.e. 
how to include the stochastic source terms in ()2.22j) . 

In this section, we shall make an explicit computation which shows that 
the form of the Green-Kubo formula, as defined e.g. in leads to the correct 
conductivity for the system with the non-zero on-site potential. 
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Theorem 5.1 Given T > 0, 

1 /'OO 1 noo 

«:gk(T) = TfrJ lim dfj\t-,T) = ^ra — dtCfj\t-T) = n. 
J- Jo J- Jo 

(5.1) 

In the theorem, k is defined by f|4.18p . and 

1 

where 



CTAt;T) = j;—j{J{qit),p{t))J{qiQlpm)s, (5.2) 



N-l N-1 n 

i=l 1=1 

The expectation value in ()5.2|) refers to the stochastic time evolution defined 
in section 121 when the initial values (qr(0),p(0)) are distributed according to 
the equilibrium Gibbs measure at temperature T. The proof is a relatively 
tedious explicit computation, which we do not report here in full detail. 

Proof: Define first the matrix K by 

{qi - q2, for z = 1 

Qn-i - Qn, for i = N , (5.3) 
Qi-i-Qi+i, otherwise 

so that 

2 ^"^ 

-^J{Q,P) = y]fe - Qi+i){Pi+Pi+i) = P^Kq. 
1=1 

Then ^^cfj\t-T) = ^g^{t)ioi 

9N{t) = y2(^^ ^ 1) WfKqit) piOfKqiO))s,. 

The initial equilibrium measure is Gaussian with zero mean and with a 
covariance C(0, 0) = TE, where E = T'^ S^"'^''^^ is by (IT^ independent of 
T. Correspondingly, 

r(, - / Te-'^E, when t > 
^^'^'UJ-I TEe'^\ whent<0 ' 
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Applying the "pairing rule" of Gaussian correlations and setting 

/C = 

we then obtain, for t > 0, 



K"" 
K 



9N{t) = ^(^^^Tr [iCe-'^ElCEe-'^''] , (5.4) 

and, for t < 0, g^it) = 5'7v(N)- 

Let us proceed by assuming the existence of lim7v-»oo gN^t) and later com- 
ment on how to prove this. From (j5.4j) we get 

l^7iv(t)|<^^^^||e-*^fl|/Cfl|i5;f. (5.5) 

Since ||$^^|| = snpj^l/fik < 1/7^; the norm of E is bounded uniformly in 
N ^ oo. The same is clearly true for ||/C||, and by ()2.1U|) and ()5.5|) we can 
now apply dominated convergence in 1)5.11) . This proves the integrability of 
the limit function, and yields 



Kqk = — lini Tr 



K / dte-'^EKEe-'^^ 



(5.6) 



We have now proved the first two equalities of the theorem. We note that 
the above argument, which allows to take the thermodynamic limit out of 
the time-integral, would fail if 7 = 0, as then neither the bound on nor 
the exponential decay of ||e~*"^|| would be uniform in N . 

Let us then denote 

S' = j^te-'^ElCEe-'''^ = ( ^^j^ ^, ) 

which is possible, as the integral clearly yields a symmetric operator. Then 
Tr [ICS'] = Tr [K^{Z'f + KZ'] = 2TT[kZ'] (5.7) 

where K = F^KF and Z' = F^Z'F. Like the matrix S defined by (fTTT]) . 
S' is the unique solution of the equation 

AS' + S'A^ = EICE=( ^"q^^ 
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In appendix^ we prove that 

^ki = 1777 r(-^-)fc«, 

where is the antisymmetric part of K, and G and Ck were defined in 
section 121 In particular, Z' is antisymmetric, and thus by (j5.7|) . 

2A ^ 1 

Tr [ICS'] = 2Tr[k_Z'] = - J] T,7-—AK-)l (5.8) 

Applying the definitions of F and i^' and neglecting all symmetric terms, 
we get after some algebra 

.f.. 2 sin(^)sin(^) 



where 5„,odd = 1, if is u is odd, and zero, if u is even. Observe then that for 
Z ~ /c we have 



{K_)ki ^ -z ^r4-«,oddSin 



7ck 



while elsewhere {K_)ki = 0{N~^). Using this observation and the equality 
Suez ^«.odd/w^ = 7r^/4, it is possible to prove that 

1 ^ r^on 8A 1 sin^(^) /I 

k=l 1^ ^bUl \^2(Ar+l) ) 

The same methods can be employed to show that limTv^oo fi'Af(^) exists 
for alH > 0. First write the trace in equation (j5.4p in the eigenspace of the 
force-matrix $, and then apply the above approximation to K^i to find that 
only terms with k I contribute and the contribution has a finite limit. 

Combining the above with equation ()5.6|) . we have now proven that 

oj"^ /"\ sin^fTTx) 

K-GK = -- / dx 



A io z/2 + 4sin2(f ) 
which shows that kqk = for all T. □ 
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6 Non-uniform heat bath couphng 



Let us now consider the case when the heat bath couphngs Aj are not all 
equal and define Ajj = 6ijXi. As long as ||e~*^|p remains integrable, we can 
repeat the computations in section |21 and conclude that equations ()2.1Hj) for 
the stationary covariance matrix are still valid. In particular, the matrix Z 
is then antisymmetric. Therefore, by redoing the computations in section 
12.11 we get the average of the energy transfer Ri and of the current Jj in the 
steady state from the equations 

{Ji)s = uj^Zi^i+i and {Ri)s = Xi (Ti - Va) . 

Thus the self-consistency condition still has the same form as before but, 
as the earlier explicit solution of the steady state covariance 5* is no longer 
possible, redoing the existence, uniqueness and local thermal equilibrium 
results is not straightforward. On physical grounds, we expect the results 
to remain valid whenever there is an A^-independent A > 0, such that the 
number of i's for which Aj > A is proportional to A^, certainly whenever this 
is true for all i. Instead of trying to redo the proofs, we shall check what 
happens if we assume that these results hold also when the Aj are not all 
equal. 

The equations for the stationary covariance (j2.13p now yield the relations 

2{Ri)s = 2cu\Zi.,,i - (6.1) 
(Aj + Aj+i)Zj_j_|_i = uj'^{Uii + f/i,j+2 — f/i-i,i+i — ?7i+i,i+i). (6.2) 

The first equation implies that the current in the self-consistent steady state 
is constant, and then, by summing (j6.2|) over i = 1, . . . , A^ — 1, we get 



„ (6.3) 



where 

1=1 

Let us next assume that the local equilibrium result proved in section ^21 
is still valid, i.e. that for every i,j 

U,,=T,<l>r^ + 0{6M), 
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where is defined by equation (j4.6j) . and that ej\f when N ^ oo. 
Choosing the Aj such that hniAr^oo ^^^^ = A > 0, we get from (|6.3|) the 
scahng of the total current, 

hm [N - 1) J(^) = k{Tl - Tr) (6.5) 

TV— ►oo 

where k is given by ()4.18|) . with A replacing A in the equation. 

The conductivity will in general be space-dependent for non-uniform cou- 
plings. Consider, for instance, any sequence of couplings X[^^ > which is 
bounded (i.e. supj ^ A^^^"* < oo) and for which the limit 

A(x) = lim — V Af ^ (6.6) 

l<i<Nx 

exists for all x G [0, 1] and defines a smooth function A with A(l) > 0. 

By summing ()6.2|) over a range of indices from 1 to j — 1 we get, after 
applying the local equilibrium assumption and ()4.22|1 . that 

2J(^) r^.(N) Af) + Al^) 



T(EM"'-^^^^)=r,-T, + oM. 

i=l 



Then, by applying ()6.3p and ()6.6|) . we can conclude that the temperature 
profile now converges to 



A(l)' 

and, therefore, that Fourier's law is satisfied with a thermal conductivity 

k(x) = = (6.7) 

A(x)2 + i/2 + ^z/2(4 + i/2) 

where X{x) = ^A(x). Note that on any interval on which A(x) = the local 
conductivity is infinite and the temperature profile remains constant. 

For instance, if Aj = A for every m:th coupling and Aj = otherwise, 
we get a finite conductivity equal to m times the one computed in section 
0] Taking m = N we then (formally) recover the linear divergence of the 
conductivity in N which was found in 0. 
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7 Higher dimensions 



Here we extend the results of the previous sections to a system of oscillators 
first in d = 2 and then also in higher dimensions. The solution can be 
obtained in a way very similar to what we did in sections El to El and we shall 
just report the necessary adjustments. Moreover, we shall only consider 
explicitly the system in two dimensions. No real modifications are necessary 
to extend the computations to higher dimensions. 
The Hamiltonian for the system is now given by 

N N' ^ 

i=i j=i 

N' N+1 N N' 

j=l 2=1 i=l j=l 

where we assume, as before, that qoj = qN+i,j = and we fix periodic 
boundary conditions in the second direction, i.e. = qi^N'- As in section 
the oscillator at the site is coupled to a Langevin heat bath with 
temperature Tij, and we set Tij = T^, Tjyj = Tr while all other Tij are to 
be determined self-consistently. Thus the time evolution is still defined by 
equations (|2.3|) - (j2.5|) if we interpret the operator A in (j2.5|) as the discrete 
Laplacian in two dimensions with mixed boundary conditions: Dirichlet in 
the first direction and periodic in the second direction. 

Extending the results in section 12.11 we first define the local energy by 

Hi,j{q,p) =^Plj + uiqij) + ^[v{qi,j - qi-ij) + v{qi+ij - qij) 

+ - Qij-i) + v{qi,j - gij+i)] (7.2) 

again with double contribution for the terms involving q^ j and gTv+ij- Then 
the analog of ()2.22j) is true if we define the current as a two dimensional 
vector with J^j = for i = 0, A^, and with the other components given by 

■Jlj{Q,P) = -yl^i+ij - Qi,j)iPi,j +Pi+i,j), (7.3) 
■JlM^P) = -yl^ij+i - Qi,j)iPij+Pi,j+i), (7.4) 
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where qij and Pij are periodic in j. 

The source terms are given by Rij = X(Ti j — (pijY), and the self- 
consistency condition thus becomes 

Tij = {iPijf)s, for 2 = 2, . . . , iV - 1 and J = 1, . . . , iV' 

with Tij = Tl and T^j = Tr. The main observation is that, as in |6j, we can 
Fourier transform this system in the periodic direction and obtain a system 
of decoupled chains. 

More precisely, let for k = 1, . . . ,N' 

^ AT' 1 ^' 

qi{k) = -j=Y,(lij(^'^''' when g^,,- = -= ^ g,(/i;)e-'^'=^ (7.5) 

j=l k=l 

and define Pi{k) analogously. This corresponds to a change to a (complex) 
eigenbasis of the periodic Laplacian, and we obtain that, for any fixed k, q{k) 
and p{k) satisfy equation ()2.Hj] with the only difference that now the potential 
$ is given by ()2.5p with u"^ replaced by i/(fc)^ = + 2(1 — cos(^)) > z/^. 
However, the noise term will then become more complicated and it can still 
a priori couple the components with different values of k. 

In general, p{k) and q{k) are complex numbers, and the stochastic equa- 
tions should be understood applying to the real and imaginary part sepa- 
rately. However, as A remains a real matrix, equation ()2.9p still holds if we 
replace the matrix by 

, where iaa\k;i',k' = 5i,2\^ T,j e'^^'^^'K (7.6) 

j 

On the other hand, our bound for the norm of the exponential of A is obvi- 
ously still valid, and we can conclude that for every temperature profile there 
is a unique stationary state which is reached exponentially fast and which is 
determined by equation 1)2.111) with the matrix ()7.(j)l replacing there. 

Next we need to prove the existence and uniqueness of the self-consistent 
temperature profile. In fact. Theorem 13.11 is valid also in the higher di- 
mensional case considered here, but since the proof remains essentially un- 
changed, we do not include it here. 

The boundary conditions we impose are constant in the periodic direction, 
and we expect from symmetry that the self-consistent temperature profile is 
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also constant in that direction, even for finite A^. This is also directly implied 
by the above quoted uniqueness since, if (tij) is a self-consistent profile, then 
also its translates, Tjj = tj.j+j' for any j', are self-consistent with the same 
boundary conditions. 

Consider thus a temperature profile Tij = ti for which ri = Tl and 
tn = Tji. By (j7.6|) . we then have always 

Applying this in the equation corresponding to ()2.12p reveals that the com- 
ponents having different values of k then become independent in the steady 
state. In particular, 

{pi(k)pj{k'y)s = {piik)pj{ky)sSkk', 
and, therefore for all i,j, 

k,k'=l k=l 

(7.7) 

Here the expectation value can be computed by {pi{k)pi{ky) s = {M{k)T)i, 
where M{k) = M\^2^^(j.-^2 and M is the matrix defined in section |21 There- 
fore, simply by replacing the matrix M with (A^')~^ Tlk-^i^) '^^^ repeat 
the computations in section El and find a vector r which leads to a self- 
consistent profile Tj J. 

It is then easy to see, applying the decoupling of the modes as above and 
then using the antisymmetry of the covariance component Z, that there is 
no average current in the second direction, i.e. {Jlj)s = 0. Similarly, we get 
for all i = 1, . . . , A^ — 1 the result 

1 

{Jl^s = u.'z,,,^i, = -J2'^'""U.ikr (7-8) 

k=l 

where J'-^^ denotes the current through the corresponding chain. 

Repeating the computations in section HI and using the above decoupling 
of the fc-modes, we can then conclude that local equilibrium holds in the 
limit N, N' oo (for this one needs to notice that the exponential decay of 
correlations is uniform in k, as i/(fc)^ > > 0), the limiting temperature 
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profile is constant in the periodic direction and connects and Tr linearly 
in the first direction. Fourier's law is also satisfied with the conductivity now 
given by 

- = ^ E = T i \ , r.... , , ...... (^-9) 



N' 



N'^ X Jo 2 + %)2 + v/%)2(4 + %)^ 



where = + 2(1 — cos(27r?/)). 

For the system with d—1 extra periodic dimensions, we could analogously 
arrive at the same conclusions, but with a conductivity 



where now ^{yY = + 2 YltZxi^ ~^ cos(27r?/j)). Observe, in particular, that 
the conductivity steadily decreases with each added dimension. We prove 
in appendix O that the asymptotic behavior of the conductivity in the limit 
d — >• cx) is given by 

where the correction term depends only on v and d. Thus by choosing a 
suitable sequence of A = 0{d~^), we can have A ^ when d ^ oo and still 
keep the conductivity finite and constant. 

It would also be straightforward to check that the Green-Kubo formula 
holds in the higher dimensional case. More precisely. Theorem 15.11 is still 
valid for the above system in d dimensions, if the current-current correlator 
is defined instead of (j5.2j) by 

Cfj\t-T) = j^y{q{t),p{t)) ■ J{q{0),p{0)))s,. 

For proving this, the ( J^)-term can be analyzed exactly as before, while the 
analysis of the remaining ( JV*)-terms in the periodic directions will be even 
simpler, as in the complex eigenbasis used here the operator corresponding 
to K will be exactly diagonal. 

8 Discussion 

We raise again the question, discussed extensively in P and [TT] , of whether 
it is possible to derive Fourier's law for a system with purely Hamiltonian 
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bulk dynamics. There are two ways of formulating this problem: (i) The 
system could be fully isolated and evolving towards equilibrium from an ini- 
tial nonuniform local equilibrium state, (u) The system could be maintained 
in a stationary non-equilibrium state by coupling it at the boundaries to in- 
finite reservoirs, either stochastically as in [Sj (or variations thereof, see ^J) 
or mechanically as in [T^ ITB] . One could also keep the temperature fixed at 
the end of the system by means of deterministic Gaussian thermostats 

In the first case this amounts to proving the existence of a hydrodynami- 
cal scaling limit on the dissipative time scale. This is a well known, extremely 
difficult problem . It is clearly not true for the harmonic crystal or other in- 
tegrable models but is believed to be true for macroscopic systems with more 
realistic type of interactions, e.g. hard spheres or with Lennard- Jones poten- 
tials. For anharmonic crystals, the kind considered in one would have 
to go beyond the Kolmogorov, Arnold, Moser domain ^Sj and presumably 
also beyond the Fermi, Pasta, Ulam models The only mechanical 
system, for which such a result has been derived, is for the highly degenerate 
model of a macroscopic system of independent particles moving in a periodic 
array of scatterers, i.e. for the multi-particle Sinai billiard, where one proves 
Fick's law, the analog of Fourier's law for the conserved particle current [T7] . 

In the second case of stationary non-equilibrium states one may hope to 
prove a global Fourier's law, i.e. we want CJc/(Tl — Tr) — > k as £ — oo. 
Here C is the distance, in microscopic units, between the boundaries of the 
system, say a cylinder, maintained at fixed temperatures and Tr. We 
want a k which depends only on the bulk properties of the system. We 
expect further that when Tl ^ Tr = T, the limit of k should coincide with 
the heat conductivity k(T) at the local equilibrium temperature T in the 
isolated time-evolving case (z). Again the only mechanical system for which 
such a result has been proven is for the degenerate system of point particles 
moving among a periodic array of scatterers where the heat current is really 
a particle current (particles pick up energy at the right wall) |T2]. The best 
that has been proven for other systems is the existence of a stationary state 
rrr.THllin] and the positivity of (T^ - Tr)Jc for fixed C [20;. 

The results proven in this paper make use of the stochastic interactions 
in the bulk to produce a local equilibrium state. This is in the spirit of the 
general work in the last two decades proving the existence of hydrodynamical 
laws in the appropriate scaling limits for systems evolving via stochastic 
dynamics [7]. We should mention here in particular the work of Kipnis, 
Marchioro and Presutti fU] who proved results similar to ours for a model 
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with purely stochastic internal dynamics. They were in fact able to consider 
a situation where the energy is strictly conserved in the bulk rather than 
just in the average as in the model considered here. This can be done also 
for a modified (more mechanical) version of their model considered by OUa 
[221 which there is an energy conserving Ornstein-Uhlenbeck type process 
producing an energy exchange between neighboring oscillators. 

The main advantage of the BRV self-consistent model is that the aver- 
age energy flow along the temperature gradient is, as seen in ()2.23|) . entirely 
Hamiltonian. As mentioned in the introduction, it might in fact be possi- 
ble to make our model entirely mechanical by coupling each site to a large 
Hamiltonian reservoir, a la Ford, Kac and Mazur [231, which would pro- 
duce an effective stochastic reservoir that would automatically, without any 
imposition of self-consistency, be at the right temperature. 

This is in fact what seems to happen effectively when we let the dimension 
of the crystal go to infinity. As shown in Appendix O after taking the limits 
t —>■ oo and ^ oo, we can let the coupling to the interior heat baths, 
which we denote by id, go to zero as d~^, and still obtain a finite value of 
the conductivity. It is clear from the analysis in Section [HI that, if we set 
Ai = Atv = ^0 and A2 = A3 = ■ ■ ■ = Aat-i = id, then the heat conductivity is 
obtained by replacing A by id in ()4.18p and in ()7.10|) . 

An open interesting problem is to consider our model for an anharmonic 
crystal, e.g. by setting in (j2.2|) . u{q) = |7^g^ + ^Sq^. We expect that for a 
fixed 6 > the heat conductivity k, would have a finite limit as the auxiliary 
couplings with the interior heat baths are taken to zero. It might even be 
possible to prove such a result by starting with a perturbation expansion in 
6 around the local equilibrium stationary state found here and then doing 
a suitable resummation or applying a renormalization group type argument. 
See however, the results of the perturbation expansion in the case with purely 
Hamiltonian bulk dynamics derived in [24 . 

We note finally that the harmonic heat conductivity k given in ()4.18p 
would remain finite if we let A — and 7 — 00 in such a way that A7^ — >■ 
a > 0. It is not clear whether this limit has any physical significance. 
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A Bound for the time-evolution matrix 

Let F be the orthogonal matrix defined by equation ()2.15|) . and define 

As F diagonalizes $, we then easily see that A = T AT^ ^ where (with Y 
again denoting the eigenvalue matrix of $) 



A 



-/ 

Y \I 



Since A is block diagonal (after a permutation of indices) and is or- 
thogonal, it follows that the norm of the exponential satisfies 

||e~*^|| = max||e-*^'=|| (A.2) 



where for each k we have defined 

Pk A 



Ak={ a (A.3) 



The eigenvalues of A^ are 



A / 
«fc = 2 ^ ^^^^^ Pk = y ^ - Pk, 



and it is easy to see that 



A 7^ 



Re ct^ > a = min| — , -r- 1 > 0. 

12 A J 



However, since Ak is not symmetric (in fact, there are values of the pa- 
rameters when it is not even diagonalizable) we have to take more care in 
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analyzing the norm of its exponential. Performing the Jordan decomposition 
of explicitly yields 

from which we straightforwardly arrive at the following bound valid for t > 0, 

||e-tA II < e-t^ [1 + t{l + 7^ + Au^ + A/2)] . 

Applying this to ()A.2j) easily yields the conclusion in sectional at equation 
fl2.10|) . We remark that if 7 = 0, we could still have a lower bound a > 0, 
but it would not be uniform in N. In fact, since then inf^ /i^ = 0{N~'^), we 
would then need to take also a = 0{N~'^). 



1 + 



iaiin.[pki) 

Du 



— Uh 



i 

-A/2 



(A.4) 



B Solution of the stationary covariance 

We derive here an explicit solution to the equation 

as + sa^ = j:^, (B.i) 

which — for the matrix used in section |2] — will yield the stationary covari- 
ance matrix. The matrix A is defined as in ()2.4j) . but we need the solution 
for a more general "noise matrix" in sectional Therefore, we consider here 

b 
2\d 

where b and d are real N x N matrices and d'^ = d. 
Denoting 

^_ f U Z 

^-{z^ V 

we get that S" is a solution to ()B.1|) if and only if its components satisfy 

Z^ = -Z 

V = -{<^U + U<^) -b+ 

2 (B.2) 
XZ = -{^U - U^) + b- 

2\{d + 6+) = $Z - Z$ + \{U<^ + <^U) 
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where 6+ and 6_ are the symmetric and the antisymmetric part of b. 
If we define 

d = F^dF and V = F^VF, 
and also b, U and Z similarly, then we get 



2\^{dki + {b+)ki) - (/Ufc - fxi)ib_)ki 
2A^(/Xfc + fii)dki - if^k - ^-i) (fJ'kbki - fJ'ikk 
{f^k - fJ'i)dki + f^khi - l^ibik 

where /i^ are the eigenvalues of and for all k and / 

Qki = Auj^G{ck, q) > 





2 


Uki = 






9kl 




1 


Vki = 






9kl 




2A 


Zki = 






9kl 



where Ck and G are defined in equations ()2.18|) and ()2.20|) . respectively. 

When 6 = and dij = 6ijTi we get ()2.1(jj) - ()2.19j) . In sectional we need to 
know Zki when d = and b = ^~^K^ with K defined by Since then 

bki = Kik/fik, where K = F'^KF, we get 



A 



1 



^kl 



uj^G{ck,ci) 

with denoting the antisymmetric part of K. 



kl 



C Asymptotic behavior of the conductivity 

In section |71 we derived a formula for the conductivity of the (i-dimensional 

crystal, 



UJ 

T 



where 



y 



[o,i]d-i 2 + i){yy + ^i^ivYiA + p{yY) 



and v{yY = z/^ + 2 X]f=i (1 ~ cos(27r?/j)). Here we prove that the asymptotic 
behavior of / for ^ oo is given for any fixed z/ > by 



^(1 + 0(1)). 



(C.l) 
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First we point out that for all r > 

1 /"^ sin^inx) /"^ sin2(27ra;) 

— da: — — r = dx 



2 + r2 + ^r2(4 + r2) Jo r^ + 4sm'^{^) Jo + 4 sin^ (tto;) 
and, therefore, 



' [0,1]'* + 4 ^^^^ sin2 (yri/i) 

Since the denominator is always strictly positive, we can then use the formula 
1/r = J^dt exp{—tr) valid for all r > and obtain 



= / dte-'''"h{t)Io{t) 
Jo 



d-l 



where 



- [ d|/e-^*^^"'("^), and 
Jo 

h{t)= I dy sin2(27r7/)e 

^0 



4f sin^ (tt y) 



Now both functions Ii{t), i = 0,1, are clearly continuous and strictly 
monotonously decreasing from 7j(0) to when t goes from to oo, with 
7o(0) = 1 and /i(0) = |. In addition, /q is bounded for alH > by 

m < (C.2) 

This follows from 

/o(i) = rdxe-^*^'°'('^^/') 
Jo 

< / dxe-"*^ / dye-"^'*/(i+*) 

~ Jo VT+t io 



since the derivative of the last integral is negative for t > 0. 
By dominated convergence we then find that, when o? — > oo, 

/oo 
dte-"''li{t)IoitY-^ 0. 
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Changing variables to s = td in the remaining integral shows then that 



Since 1 + x > 2"^ for all < x < 1, inequality ()C.2|) yields the bound 



for all < s < c? and d > 2. This means that dominated convergence can 
also be applied to the integral in ()C.3|) . and as Io{s/d) = 1 — 2s/d + O^d'"^), 
we then find 



which proves the equation (jdlll . 

References 

[1] F. Bonetto, J. L. Lebowitz, and L. Rey-Bellet, Fourier's law: a challenge 
to theorists. In A. Fokas, A. Grigoryan, T. Kibble, and B. Zegarlinski 
(eds.), Mathematical Physics 2000, pp. 128-150, London, 2000. Imperial 
College Press. 

[2] S. Lepri, R. Livi, and A. Politi, Thermal conduction in classical low- 
dimensional lattices, Phys. Rep. 377 (2003) 1-80. 

[3] M. Bolsterli, M. Rich, and W. M. Visscher, Simulation of nonharmonic 
interactions in a crystal by self- consistent reservoirs, Phys. Rev. A 4 
(1970) 1086-1088. 

[4] M. Rich and W. M. Visscher, Disordered harmonic chain with self- 
consistent reservoirs, Phys. Rev. B 11 (1975) 2164-2170. 

[5] Z. Rieder, J. L. Lebowitz, and E. Lieb, Properties of a harmonic crystal 
in a stationary nonequilibrium state, J. Math. Phys. 8 (1967) 1073-1078. 

[6] H. Nakazawa, On the lattice thermal conduction, Suppl. Progr. Theor. 




(C.3) 





Phys. 45 (1970) 231-262. 



31 



[7] H. Spohn, Large Scale Dynamics of Interacting Particles. Springer, 
Berlin, 1991. 

[8] B. 0ksendal, Stochastic differential equations: an introduction with ap- 
plications. Springer, Berlin, fifth edition, 1998. 

[9] J. L. Lebowitz and H. Spohn, A Gallavotti- Cohen type symmetry in 
the large deviation functional for stochastic dynamics, J. Stat. Phys. 95 
(1999) 333-365. 

[10] G. L. Eyink, J. L. Lebowitz, and H. Spohn, Hydrodynamics and fluc- 
tuations outside of local equilibrium: Driven diffusive systems, J. Stat. 
Phys. 83 (1996) 385-472. 

[11] J. L. Lebowitz, Exact results in nonequilibrium statistical mechanics: 
Where do we stand?, Suppl. Progr. Theor. Phys. 64 (1979) 35-49. 

[12] J. L. Lebowitz and H. Spohn, Transport properties of the Lorentz gas: 
Fourier's law, J. Stat. Phys. 19 (1978) 633-654. 

[13] J.-P. Eckmann, C.-A. Pillet, and L. Rey-Bellet, Non-equilibrium sta- 
tistical mechanics of anharmonic chains coupled to two heat baths at 
different temperatures, Commun. Math. Phys. 201 (1999) 657-697. 

[14] G. Gallavotti and E. G. D. Cohen, Dynamical ensembles in stationary 
states, J. Stat. Phys. 80 (1995) 931-970. 

[15] M. Tabor, Chaos and Integrability in Nonlinear Dynamics: An Intro- 
duction. Wiley, New York, 1989. 

[16] E. Fermi, J. Pasta, and S. Ulam, Studies in nonlinear problems, I. In 
A. C. Newell (editor). Nonlinear Wave Motion, pp. 143-156. American 
Mathematical Society, Providence, RI, 1974. Originally published as Los 
Alamos Report LA-1940 in 1955. 

[17] J. L. Lebowitz and H. Spohn, Microscopic basis for Pick's law of self- 
diffusion, J. Stat. Phys. 28 (1982) 539-556. 

[18] S. Goldstein, J. L. Lebowitz, and E. Presutti, Stationary states for a 
mechanical system with stochastic boundaries. In J. Fritz, J. L. Lebowitz, 

and D. Szasz (eds.), Random Pields (Colloquia Mathematicae Societatis 
Jdnos Bolyai 27), pp. 403-419, Amsterdam, 1981. North-Holland. 



32 



[19] S. Goldstein, C. Kipnis, and N. laniro, Stationary states for a system 
with stochastic boundary conditions , J. Stat. Phys. 41 (1985) 915-939. 

[20] L. Rey-Bellet and L. E. Thomas, Fluctuations of the entropy production 
in anharmonic chains, Ann. H. Poinc. 3 (2002) 483-502. 

[21] C. Kipnis, C. Marchioro, and E. Presutti, Heat flow in an exactly solvable 
model, J. Stat. Phys. 27 (1982) 65-74. 

[22] S. 011a. Private communication. 

[23] G. W. Ford, M. Kac, and P. Mazur, Statistical mechanics of assemblies 
of coupled oscillators, J. Math. Phys. 6 (1965) 504-515. 

[24] R. Lefevere and A. Schenkel, Perturbative analysis of anhar- 
monic chains of oscillators out of equilibrium, preprint (2003), 
fhttp : // arxiv . org/abs/math-ph/0303050| 



33 



